% This file is used to define the Global variables, initial conditions and call the required files.
% The variables are:
% p: represents zeta : wave height
% q : represents phi : surface velocity potential
% 1-D case .. in X direction


%--- main parameters
global U  h  g  p  q  nt dx h1 x tmax L steps y itermax Q p_global Length;
range=[ 0.5] ;
for k = range
    %clear all;
    
    %----- Physical Parameters----------
    
    U = 1; % Current Velocity of the wave. Assumed Constant
    h = 50; % height of the bed, assumed constant
    g=  9.81; % Gravitational Constant
    
    %index = find(range==k);
    %------------ Simulation Parameters % -----------------
    tmax = 100;     % max time, s
    Length=  100;
    x = 200;         % number of axial nodes
    steps = x-1;    % number of axial steps
    dx = Length/steps;   % axial step spacing, delta x
    % Check the CFL condition:
    dt1 = 0.5*dx/sqrt(g*h);
    %--- Initial Conditions ----------
    pi1=1;
    qi=1;
    
    for i=1:x
        p(i,1) = pi1*cos(4*pi*(i-1)*dx/Length);
%        p(i,1) = 0.25*pi1*cos(4*pi*(i-1)*dx/L) + 0.25*pi1 *cos(8*pi*(i-1)*dx/L)  + 0.25*pi1*cos(16*pi*(i-1)*dx/L)+0.25*pi1*cos(32*pi*(i-1)*dx/L)  ;
        q(i,1) = qi;
    end
    
    % Run the matlab based Codes
    i=1;
    
    Complete_Matlab(2);
    A(1:x,i)=y(end,1:x);
    B(1:x,i)=y(end,x+1:2*x);
    % i=2;
    % Complete_Matlab(i);
    % A(1:x,i)=y(end,1:x);
    % B(1:x,i)=y(end,x+1:2*x);
    
    % Run the custom Program
    %dt =0.023;
    %CFL = [14.4 7.2 3.6 1.8 0.9 0.45];
    %CFL=1.9;
   % dt_v = CFL*0.5*dx/sqrt(g*h);
    %dt_v=[CFL*dt1];
    %dt_v = dx*50;
    dt_v = [0.1];
    for i = 1:length(dt_v)
        dt = dt_v(i);
     %   X = sprintf(' Simulation with Nodes: %f , CFL %f ,  dt %f',x,CFL(i),dt);
     %   disp(X);
       TimeDiscretizers(dt,1);
        A(1:x,i+1) = p;
        B(1:x,i+1)=q;
        X = sprintf('The Norm of the Error for Zeta is %f' , norm_L2(A(:,1),A(:,i+1)));
        disp(X);
        X = sprintf('The Norm of the Error for varphi is %f' ,  norm_L2(B(:,1),B(:,i+1)));
        disp(X);
    end
   %     A_G(:,:,index) = A;    
end
